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An integrable case of the p + ip pairing Hamiltonian 
interacting with its environment 

Inna Lukyanenko, Phillip S. Isaac, Jon Links 


Abstract 

We consider a generalisation of the p + ip pairing Hamiltonian with external interaction terms. These terms 
allow for the exchange of particles between the system and its environment. As a result the u(l) symmetry 
associated with conservation of particle number, present in the p + ip Hamiltonian, is broken. Nonetheless the 
generalised model is integrable. We establish integrability using the Boundary Quantum Inverse Scattering 
Method, with one of the reflection matrices chosen to be non-diagonal. We also derive the corresponding Bethe 
Ansatz Equations, the roots of which parametrise the exact solution for the energy spectrum. 
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1 Introduction 

Understanding how quantum systems interact with their environment, and being able to control such interactions, 
is a major challenge facing quantum engineering. One such framework where this applies is provided by Josephson 
junctions, fabricated through weakly-coupled superconductors. These structures have received widespread study as 
a potential architecture for the coherent control of quantum bits e.g. [1-5]. In some instances, such as the Cooper- 
pair box Josephson junction, the system is described in terms of a simple Bose-Hubbard tunneling model [1,2]. 
More refined analyses to produce insights into environment interactions have also been undertaken, in particular 
through explicit use of the degrees of freedom of the s-wave pairing Hamiltonian in the strong-coupling limit [6,7]. 

In recent years the p + ip pairing Hamiltonian has emerged as an example of a superconducting model which 
is integrable, and admits an exact Bethe Ansatz solution [8]. The solution was obtained by application of the 
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Quantum Inverse Scattering Method (QISM) [9,10] associated with the trigonometric XXZ solution of the Yang- 
Baxter equation in the quasi-classical limit. This result has led to several studies of the model, including re¬ 
derivations of the solution from different perspectives such as the Richardson-Gaudin approach, investigations into 
the ground-state structure, and extending the application of the exact solution for the calculation of correlation 
functions [11-17]. 

Here we will establish that there is an extension of the p+ip pairing model, involving interaction terms coupling 
to the environment, which maintains integrability. These interaction terms allow for the exchange of particles 
between the system and its environment. As such they break the u(l) invariance associated with conservation of 
particle number which is present in the p + ip Hamiltonian. From a naive perspective it appears that the prospect 
for constructing an integrable extension of the p + ip pairing Hamiltonian with broken u(l) symmetry is dire. The 
p + ip pairing Hamiltonian is constructed as a linear combination of mutually conserved operators with co-efficients 
which are dependent on the particle number [11-13,15]. If the particle number is not conserved, this approach 
fails. Fortuitously, our recent study [18] uncovered a means to overcome this issue. It turns out that the p + ip 
pairing Hamiltonian and its exact solution can also be derived through application of the Boundary Quantum 
Inverse Scattering Method (BQISM) [19] associated with the rational XXX solution of the Yang-Baxter equation 
in the quasi-classical limit. In this approach, the Hamiltonian is constructed as a linear combination of mutually 
conserved operators with co-efficients which are independent of the particle number. Thus this framework can be 
extended to produce a generalised model interacting with the environment, as we will describe below. 

In Sect. 2 we introduce the Hamiltonian and review the Boundary Quantum Inverse Scattering Method of 
Sklyanin [19]. Sect. 3 details a construction to obtain a set of mutually commuting operators, and shows how the 
Hamiltonian may be expressed in terms of the elements of this set. Sect 4. then utilises Bethe Ansatz results to 
obtain the exact solution for the Hamiltonian. Concluding remarks are given in Sect. 5 


2 Preliminaries 

2.1 The pairing model interacting with its environment 

We first introduce the isolated pairing mode not interacting with the environment. Let Ck, cj^. denote the anni¬ 
hilation and creation operators for two-dimensional fermions of mass m and momentum k = (kx,ky). Then the 
pairing Hamiltonian is 


Hn = 


k 


[1^ 

2m 


cj^Ck 


G 


(kx + iky)(k'^-iky)clclj^c-k^ck> 

k^±k' 


where G € M is a constant and the summation is taken over all momentum states k. The annihilation and creation 
operators Ck,cj^ satisfy the canonical anticommutation relations: 

{Ck,Ck'} = {Ck’4'} = {Ck.Ck/} = 4k'1- 


Now consider a more general Hamiltonian with an extra term 


H = Hn 


2 ^ 
k 


(^(4 -I- iky)clcl^^ + (4 - ffcy)c_kCk) , 


( 1 ) 


where F € M is a constant. We note that this Hamiltonian is Hermitian, and the extra term can be interpreted as 
creation and annihilation of pairs of fermions, resulting from interaction with the environment. It is important to 
distinguish this type of interaction with the environment from other examples, e.g. [20] in the context of a heat 
bath, which facilitate a notion of entanglement with the environment. In our model there is no entanglement 
between the system and the environment, because the state space for the environment is not explicitly defined. 
We will comment further on this aspect in the Conclusion. 

We now restrict to the Hilbert subspace that allows only paired particle states. By imposing this restriction, 
we do not consider states on which the operators in the interaction term in the Hamiltonian (i.e. the second term) 
has trivial action. On this subspace the following equality is satisfied: 

2cj^CkcL^c_k = clck + cLi^c_k. (2) 
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Set Zk = |k| and k^; + iky = |k|exp(z0k)- Introduce the following notation: 

= exp(z0k)CkcL^, = exp(-i0k)c_kCk, = cj^cLj^c.kCk - 
Remark 2.1. On this restricted subspace, one may verify the su(2) algebra commutation relations: 

[S+,S^] = 2S^. 

We now use integers to enumerate the unblocked pairs of momentum states (k and —k). Working in units such 
that m = 1, using equation (2) and ignoring the constant term - we obtain 

fc=i 

c c 

Ho = Y. - GE E 

k—1 k—lj^k 


which exhibits u(l)-symmetry associated with the operator = E! ^k- The full Hamiltonian that we work with 

fc=i 

therefore becomes 

c 

H = Ho + TY,Zk{S+ + S^). (3) 


fc=i 


This Hamiltonian no longer possesses u(l)-symmetry. 

In this article we show that the Hamiltonian (3) is integrable by means of the BQISM. Recently, a systematic 
method, referred to as the Off-Diagonal Bethe Ansatz (ODBA) has been proposed for solving such models [21]. 
This method has been since applied to several long-standing problems and the results has been summarised in the 
recent book by Wang et al. [23]. Based on the results from [22], we derive the formulae for the eigenvalues of the 
conserved operators, the corresponding Bethe Ansatz Equations (BAE) and the energy (i.e. the eigenvalue of the 
Hamiltonian). 


2.2 Boundary Quantum Inverse Scattering Method 

In this section we review the Sklyanin’s BQISM [19] and specify the ingredients in the context of our model. 
Throughout this paper we fix a vector space V = C"^. A key element of the BQISM is the i?-matrix, which is an 
invertible operator R{u) G End(y 0 V) depending on a spectral parameter u G C and satisfying the Yang-Baxter 
Eguation (YBE) in End(E 0 E 0 E) 


Ri 2 {u - w)i?i3(M)i?23(^^) = R23{lj)Ri3{u)Ri2{u - v). 

Here, as usual, the subscripts indicate the spaces in which the corresponding i?-matrix acts non-trivially. 

In this paper we consider the rational i?-matrix that is usually associated with the XXX spin chain 

-G 7] 0 0 

u rj 0 

rj u 0 ’ 

0 0 u-GrjJ 

where rj G C is the quasi-classical parameter and P is the permutation operator in E 0 E. 

In the BQISM framework the boundary conditions are encoded in the left and right reflection matrices, or 
K-matrices, K~{u) and K'^{u) G End(E), which satisfy the reflection equations in End(E ® V) 

Ri 2 {u — v)Kf {u)R 2 i{u -G v)Kf (ri) = Kf {v)Ri 2 {u + v)Kf {u)R 2 i{u — ri), (4) 

Ri 2 {v - u)K^{ u)R 2 i{-u -V- 2r])K:^{v) = K:^{v)Ri 2 {-u - v - 2ri)K^{u)R 2 i{v - u). (5) 


R{u) = ul -G rjP = 
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One can check that the following i^T-matrix satisfies the first reflection equation (4): 


K (u) 


/C + U Ip u \ 

V 4>~u C" - uj 


Then, 

K (u) = K ( U "/y) — >■!/;+,0-|—>-0+ 

automatically satisfies the dual reflection equation (5). 

We may express the Hilbert space of states in the form 


f^^+u + ri ilj+{u + r])\ 
y (/)+(« + 77 ) ^+ —u —ryy 


c 

n = (^Vj = v‘^^, 
i=i 


( 6 ) 


where each local space Vj (a copy of V) is a fixed representation space for the su(2) algebra spanned by S~ ,S^ 
(indices indicate in which space the corresponding operator acts non-trivially). For each label j in the tensor 
product (6), we introduce the Lax operator 


Laj{u) = 


I (u + rjS^ r]S, 


r]S. 


u — pS", 


= / + 


St -s\ 


SJ 


e End(K ®yj), 


where the auxiliary space Va is another copy of V. 

It is straightforward to check that the Lax operator (7) satisfies the RLL relation in End(Va ®Vb®Vj) 




(7) 


( 8 ) 


where 14 = E is another auxiliary space. 

Remark 2.2. Note that the Lax operator (7) satisfies the following property: 

Laj{u)Lajiv - u) = 

where sj is the value of the spin on the local space Vj. 

Define the monodromy matrix as 


7a(^) ^c') • • • ^l); (70) 

where Sj G C are the inhomogeneity parameters. Erom the RLL relation (8) it follows that the monodromy matrix 
(10) satisfies the RTT relation in End(Ea ®Vb®T-L) 


Rab{u - v)Ta{u)Tb{v) = Tb{v)Ta{u)Rab{u - v). 


( 11 ) 


Let us construct the dual monodromy matrix as 

Tain) = Lai{u + ei + p). ■ ■ Lac{u + ec +v)- 

From the property (9) of the Lax operator it follows that 

faiu) (xL-t{-u-ei)...L-^{-u-ec) =T-\-u), 

which implies that Ta{u) satisfies the following relations: 

Tb{v)Rabiu + v)Taiu) = Ta(u)Rab{u + v)fb{v), 
fa{u)fb{v)Rab{v -u) = Rab{v - u)fb{v)fa{u). 
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( 12 ) 

(13) 



Now, the double row monodromy matrix is constructed as follows: 


Ta(u) = Ta{u)K^ (u)fa{u), 

and the relations (11), (12) and (13) imply that it satisfies 

Rab{u - v)Ta{u)Rba{u + v)Tb{v) = Tb{v)Rba{u + v)Ta{u)Rab{u - v). (14) 


The double row transfer matrix is then defined as 

t{u) = tro {K'^{u)Ta{u)) . (15) 

Using (14), one can show that these transfer matrices (15) commute for any two values of the spectral parameter: 


[t(u), t(u)] = 0 for all u,vGC. 


Thus, (15) can be used as a generating function for the conserved operators of the system. 

In what follows, it is convenient to make a variable change u ^ u — rj/2, Sj ^ ej—ri/2 and redefine all functions 
taking this into account. This results in 


K (u) 

K^{u) 


+u-r]/2 

\(j)-{u-r]/2) 

+U + 77/2 

\(j)+{u + T]/2) 


{u-ri/2)\ 
-u + ril2) ’ 

7/>+(-« + ? 7 / 2 )\ 

-u-r]/2) ’ 


Ta{u) = Lal{u + Si) . . . Lac{u + ec)- 


Thus, the transfer matrix (15) will take the following form: 


(16) 

(17) 


t{u) = tia {K+{u)Lac{u - Ec) ... Lai{u - ei)K^ {u)Lai{u + ei)...Lac{u + ec)) ■ (18) 


3 Construction of the conserved operators and Hamiltonian 

As discussed in the introduction, our focus will be on taking the quasi-classical limit 77 —)■ 0, thus connecting our 
study to the Richardson-Gaudin class of models [13,15,16,22,24,25]. Indeed, the resulting model we refer to as 
the open, rational Richardson-Gaudin model. To be able to take this quasi-classical limit, however, we require that 
the AT-matrices satisfy the following condition: 

K~^{u)K~{u) ^ f{u)I as 77 —>• 0. (19) 


Assume the following dependence of the parameters on rj: 


C = C + Ip+ = if (/)+ = (/) -p rj\, 

= -(,+ 77 / 3 , 'ip~ = Ip + rjd, (j)~ = (j) + rjpL. 

Now consider 

^ (f;; „). ,.^(1+- ew 

Thus, the condition (19) is satisfied. Now, expanding the AT-matrices in 77 we obtain 

A'+(u) = {u) -\- pK^iu) o{ri) 


with 


K+{u) 


Pf-i-U ^fu \ 
\ (j)U f-uj 


K+ {u) 


f a -h 1/2 ju -h 7 /i/2 \ 

[Au-i-p/2 a-1/2 J ’ 


( 20 ) 


( 21 ) 
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and 

with 

Ki {u) = 

For the Lax operator we have 


K (u) = (u) + rjK^ (u) + 0 ( 77 ) 

K 2 (u) = 


( 22 ) 


+ u tpu 

(j)U — u 


13-1/2 6u-'tp/2\ 

fxu-(l)/2 /? + 1/2 j ■ 


TJ 

— 3 ^aj^ with £(ij — ( 


s: 


SJ -S^ 


(23) 


3.1 The first family of conserved operators 

In the quasi-classical limit, the conserved operators tj are constructed as follows from the transfer matrix (18): 

lim {u — ej)t{u) = rf'Tj + o{ri^). 


Substituting (21), (22) and (23) into (18) we obtain 
lim {u — £j)t(u) = 77 tr, 


^ta )^a3 ^la i^j) + V 

k>j 


^la {Sj)^ak^aj i^j ) 

^7 


k<J 

+ ’lf: + vKt{^,}e.,Ki.{e,) 

k=l 


3 ( 77 ^). 


One can check that 

and 

Thus, we have 


{Klai^j)^ak^ajKla{^j)) ~ (-^la )f (^7 )) ■ 

(T* - ^ I ^ 




£j Sk 


fc=l 


£j + £k 


+ tra {K+^iSj^ajK^^iej)) + tr^ {K+^{ej)eajK^aisj)) ■ 


Compute the traces: 


tr, {Klie,)£akia,K^aiej)) = ((1 + - C") {2S^S^, + S+S^ + SJS +), 

tr, (i^+(£,)4,1 ^i-(£,) 4fc) = 2(8, + Oie, - e)5|5/ - (e, - O^S+S^ - (sj + 0"S-S+ + 

+ 23p£j ({sj + ^)Sj S// + {sj — S/.) + 2(j)£j ({sj + ^)Sj SI + {sj — ^)S/ S^.) -\- 

+ + <t>^S-SZ - 2i3(j>S/Sl) , 

tra {Kta(.^j)£ajK:[^{.ej)) = {2a£j - C + (At/ - l4>)e]) S'/ + (^{ai/ - -f^Sj - + jsfj S++ 

+ {{a/) - XC)£j - - Aef^ SJ, 

tr, {KtAe3)£a,K^a{ej)) = {2pe, - C + (# - i3p3)e]) + (^(/3t/ + - fc?) S+ + 

+ ({134’ + ^k^)8j ~ ■ 
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The sum of these four terms leads to a family of conserved operators for the open, rational Richardson-Gaudin 
model: 

(1 + ~ r 4- 4-n 

=11 . _ * + s/s* + SJS/] + 


k^j 

C 




+E^— 

£j ^ Cf- 




2{e, + Oisj - OS]SI - {e, - O^SfS^ - {e, + O^S-S++ 


(24) 


+ {{sj + OS]S^ + (sj — OSj'S]) + 2(j)ej ({sj + OSj S] + {sj — OS]Sj .) + 

+ {rs]-s+ + O^S-SZ - 2^<PS]S]) 

+ [2{a + j3)£j -2^ + ii){\ - Os] - 4>{l - 0^“]] S]+ 

+ [tp{a + P)ej - 0l - 0^3 - + (7 - Os'j] S]~ + 

+ [(j){a + j3)£j - ^(A - fj,Oj - - (A - /r)e|] S~. 

3.2 The second family of conserved operators 

Note that we have only considered one of two possible families of the conserved operators. The second family is 
constructed as follows from the transfer matrix (18): 

lim {u + ej)t{u) = + o{ri'^). 

U^-Cj 

Here we show that these are equivalent conserved operators, i.e., fj = —tj. Let 

t{u, i) = tVa (^K]; {u)Lac{u - £ c) ■ ■ ■ Lal{u - £i)K~ {u)Lal{u + £ 1 ) . . . Lac{u + S/;)^ . 

Consider t{u, 0'^, where T = ti .. .tc denotes a transpose over all spaces. Using 

(tr,H,)‘-‘^ = tr, = tr, , 

the fact that Lax operators are symmetric 


Lajiu) — 


fu + T]{S]V 

vis+V \ 

1 /u + rjSj rjSj A 

1 ^(^7)^ 

u-v{S]VJ 

u \ rjSj' u - rjS]) 


— Laj(u), 


and an observation that K'^{u)'^ = K^{u)\^+^ff,+ and K (m)^ = K (u)|^-<_).0-, we obtain 

= iia[Lac{u ^ £c) • ■ • Lai{u ^ ei)K~{uf Lai{u - ei)... Lac{u - £n)K:^{ u)^^ = 

“ 1 “ £C ) * • ■ -^al ('^ “ 1 “ £l)^a ('^ ) ■ ■ • ^aC ,ip~ ■<r^(p~ 

= t{u, 0\lp+<£nf>'^ , 1 p~<£nf>~ ■ 

Thus, we obtain the following equality: 

t(u, 0 — 1(^5 0 . 

It follows that 


Thus, we have 


lim {u + £j)t{u,£) = lim [u - {-£j))t{u,-£) ■ 




Computing fj(£) = from (24) we obtain that fj(0 = —Tj(£)- 
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3.3 The case when one h'-matrix is diagonal 

Hereafter, we will only consider the spin-1/2 representation of this algebra acting on V: 


S+ = 





Remark 3.1. Note that in this instance Laj{u) = u ^Raj{u — rj/2). 

It now turns out that six of the parameters appearing in (24) are superfluous and can be eliminated by 
appropriate basis transformations and redefinitions of variables. First note that we can set j3 = 0 without loss of 
generality, since the dependence of (24) on a and j3 is only through the sum a j5. Next, the Lax operator is 
invariant under the local basis transformations, i.e. 

XaX,Laj{u)X-^X-^ = Laj{u) 

for any invertible X e End(C^). Thus we can almost always choose a basis in which one of the iF-matrices is 
diagonal. (The case when a iF-matrix is not diagonalisable has been discussed in [24]). For our purposes, we 
assume that K~{u) is diagonal, so that 


K (u) 
K^{u) 


+u-ri/2 

[ 0 

+ M + ri/2 
\(j)+{u + T]/2) 


C -u + tj/2)^ 
iI;+{u + t]/ 2)\ 
-u-'nl2j ’ 


For the expansion (20) this means that ^/’ = (/) = (5 = /x = 0. Substituting these into (24) we obtain 


'i — ~ + 0 


c 

E 


E 


1 


1 


Sj — Ek Ej + Ek 
1 £j 


2SJSI + J2 


1 Ej-^ 


Sj — Ek ej-\-ekSj-\-^) ^ ^ 


StS7 + 


Ej — Ek Sj H- Ek Sj 


2aEj 


- 


2 ^ 




;S! 


—() 7 7 - 7-—E 7 - et—E + 

7£i 


^ 73 - i 


Finally we may set ^ = 0 without loss of generality, although this is more technical to establish. Using the 
properties of the spin-1/2 representation, namely 

S+S- = \l + S-S+ = \l- S\ = E, 


and the identities 


-3 ^3 

we may simplify the expression for Tj to write 


1 1 Ej — ^ _ 2Ej{Ek+0 

Ej-Ek Ej+Ek£j+^ iEj-El){Ej+0’ 

1 1 Ej + ^ _ 2Ej{Ek — ^) 

Ej-Ek £j+£k£3-^ {e"^, - El){Ej - ^) 


SjTj 






-rN^sjsi+ j:^ 


2 ^j fsk+^„j.„- £k ^ o-( 


k^j ''k 

2 a£^ 


7,+i 


SfSk + ) + 

3 k 3 ' 


£j-i^ ^ £j + £, 


3 ^ ^ ^3 

Ej-S, ^ 4 


2E^j-e- 


(25) 




Consider the following local transformation on the jth space in the tensor product: 


Uj = diag 


? 

Sj-C 


Under these transformations we have 


U,S]U-^ = S], 


u^s^u-^ = 


^ q+ 


u,s-u-^ = 


o- 

ej+C^' 


Define 


.( 1 ) 


4£" 


2 £? 


k^j S' 


k^j s - ^k 


^k + ^ _ 


St S7 + ^ C'- 0 + 


Ej-^ 


2ae] 

sjst 1 + 


e]-e ^ 


^2 

I /'C. A 

+ -^.St - ^—S-. 


7 £| , Ae? 


-r C.J 

We see that, up to a constant term, it is the same expression as (25). Under the global transformation U = 
U 1 U 2 ■ ■ - Uc we define 


rf ^ = C/rj^^U-i = 

" 45? 




2aej 


= E + S-St) + -^^s: 


k^j S ^k 


+ 


7e- 


£| - e 


E- 


k^j s s 

Ae^ 

^=^J- 


£ 


S - 


2 3 


Next simply rescale to obtain 


(3) ( 2 ) 


S = 


2 S = 


^ 4fe2-£2) £ ‘^\l7‘]-e^/£l-e 

= E [2^2 S!SI + ^ - [StS- + S-St) + 2a5|+ 

k=^j 


p2 _ 2 

i fe 


p2 _ 2 


^js'i-esf-xjsj-es. 


Now we apply a change of variables Sj !->■ to obtain 


7 =E + E + 7 'N) + 207+ 


k^j S 


+ 7£jS'+- A£j5' - = 


£jTj 


i^j — 0(7 j + 0 


-I. 


5=0 


(26) 


This affirms that we may also set 5 = 0 without loss of generality. 

We refer to the set of mutually commuting conserved operators {r* : j = !,...,£} as the open, rational 
Richardson-Gaudin system in the spin-1/2 case. Note that the coefficients of the SjS^ terms in (26) are not 
antisymmetric with respect to the interchange of indices j and k. This distinguishes this set of commuting 
operators from those obtained by the Gaudin algebra approach [13,15,16] 
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3.4 Hamiltonian 

Let us now construct the Hamiltonian (3) from these conserved operators. Consider 


i: +«-»+) +2a 


j,k:j<k 

i=i i=i 


Setting A = — 7 , and making the change of variable Zj = £ we obtain 


1 


— V£-V = 

2a^ ^ C 
i=i 


j = l j,k:j<k j = l 

We see that (27) is equivalent to (3) by identifying a = G~^ and 7 = 2rG“^. 

4 Eigenvalues, Bethe Ansatz Equations and the energy spectrum 

We now turn to investigating the eigenvalues of the conserved operators, making use of the results of Wang et 
al. [ 22 ]). 


4.1 Eigenvalues 

Rewrite the il-matrices (16) and (17) in the following form (using the notation from [22]): 

K~{u) =C+ iu-v/‘2)hi-a with hi = 2 ^ ’ 


K~^ (u) + {u + r]/ 2 )h 2 ■ (f with ^2 = 


+ (j)+ i('0+ - (j)+) 


where a = (tr^, ,a^) are the Pauli matrices. To match the notation from [23] we need to normalise the vectors 

hi and ^2 

P = P = 

^ x/T^-r + i’ ^ 

The RT-matrices can then be written as 


K {u)=^/'ip (j) + 1 I ^ ^ + (u- 7 ?/ 2 )h)’■ a 


K^{u) ='\/^+(/)+ + 1 


^'tj)+(j)+ + 1 


ri/2)hl 


Let {vk I fc = 1 , 2 ,... £} denote a set of parameters that will be utilised to determine the eigenvalues of the transfer 
matrix (18). From [22], the formula for the eigenvalues of (18) is 

A(m) = + ly/'tlj+(j}+ + 1 q,(^) Q(^ + ^) h) _^^(^_.^^2)(u + ?7/2)^^ , 

Q[u) Q[u) Q{u) 
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where 


Qiu) = - Vi){u + V,), 


i{u) = 


d{u) = 


2u — rj 
2u 

2u + r] 


r 


2u 


u — 


\/ip (p +1 

r 

y/lp-(j)- + 1 


+ W 2 




\/v>+(/)+ +1 


+W 2 n 


- 77/2 u - 




Z = 1 

c 


y/'lp+(l)+ + 1 


-^/2 n 


Z =1 


{u-ei- v/2){u + ei- r]/2) 
{u-si){u + ei) 

{u- £1 + r]/2){u + ei+ 7?/2) 
{u-£i){u + £i) 


F{u)= n 




(m + - r]/2){u - £i- r]/2){u + £i+ 7 ?/ 2 )(m - £» + 77 / 2 ) 

(M-£i)(M + £i) 


c = 2 


(s; 


0 ro 




The constant c can be computed as 


c = 2 


hi ■ 


^/ip-cp- + i^/i^+(p+ +1 


-1=2 


h{'>p (p+ +(p + 1 

yJ-lp-Cp- + lx/^+<p+ + 1 


- 1 


Finally, we obtain 

A{u) =\/ + 1 


c c 

, lYT iu-Vi + v){u + v^+T]) ..T-r {u-Vi-r]){u + Vi-'n) 


+ c{u^ - r]^/4:) 


((■u + £j )2 - ?72/4)((m - £i)2 - 77 ^/ 4 ) 
( 772 _£ 2 )(^ 2 _^ 2 ) 


where 


i(M) = 


(i(M) = 


2u — rj 
2u 

2u + ri 
2u 


r 


\/^ 4> +1 2 


— 


r 


u — 


+^1 fr 

i/pp+pp+Ti 2y 

t+ ^ 




^Tp-cp- + 1 2y \/'ip+(p+ + 1 2y (u-£;)(« + £;) 


n 


{u-£i- T]/ 2 ){u + £i- 77/2) 

(7i-£i)('U + £i) 

{ u -£ i + 77/2) (m + £; + 77/2) 


(28) 


4.2 Quasi-classical limit of the eigenvalues 

The eigenvalues in the quasi-classical limit (77 —)• 0) are constructed as follows: 

lim {u — £j)A{u) = rf\j + o{rf'). 

u—>ej 

We compute this limit assuming the same dependencies (20) as for the conserved operators: 


c 


lim {u 

u->-ej 


foOlI 

1=1 


{u + £i- 77 / 2 ) (m - £; - 77 / 2 ) 
{u - £i){u + £ 1 ) 


C 


lim (u 

u=,ej 


foOlI 

1=1 


{u + £i+ 77 / 2 ) (m -£i+ 77 / 2 ) 
{u-£i)(u + £i) 







i/ip-cp- + 1 =s/'<p(p+l + oir]), 

y/lp+(p+ + 1 =y/lP(P+l ^ 1 -P + 0 ( 77 ). 
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Consider the expansion of the first two terms in (28) up to the second order in rj: 


2/^2 

+ ly/^+(l)+ + la{u) = 1 ) ~\{^ + 

+ y((A + Ai)V' + (7 + < 5 )«^') + ia + P)ej^/'^|J(|)+l - ^-^|^=((A - /x)')/'+ (7 - - ?(« 


-P) + 


+ \{e^,w+i)-e)(:^+i2( 

V ’ k^J ^ 


£j — £k £j + £k , 


0{Tf), 


2 / 

+ lyj'4)+(l)+ + ld{u) = |(£f(^<))+ 1 ) + y ( “ ^ “ yW + 1 ) + 


+ ^((•^ + ^k)i’ + (7 + 5 )(t>) - {a 


+ \m^+^)-e) (^+e( 

\ k^j ^ 




2‘\/'ip(j) + 1 


£j ~ £k £j + £k, 


((A - + (7 - 5)(j)) -?(«-/?) + 


3 ( 772 ). 


Also 


ft ‘‘‘ 7'"+"ii" 1 +y (^+^ 

J--*- {U-Vi){u + Vi) ^,\U-Vi U + Vi 

t—1 ^ ^ ' 1—1 ^ 


{u - Vi - ri){u + Vi - ri) ^ 

(1! - IT ■ 1 C/ 4- 71A ^ f ^ 


{u-Vt){u + Vi) 

Combining these calculations then leads to 


1 


1 


U — Vi U + Vi 


0(7), 

0(77). 


lim {u-£j)^Jil} (j) +l^/tp+(j)++ la{u)Y[ ^^ + r]){u+ v^ +-q) ^ 

[U — Vi)[UVi) > 


i^l 


77 

y 


2 

c 


- {£]{'iP4) +1) - e) E yE + +1) - E yE’ 

7 = 1 ^ ki^J^j ^k 


- ^(V'((' + 1) - ^ - if((A + 7^)^^ + (7 + 'J)<(’) - {a + P)£jy/^p(t)+l+ 


i£j 


‘2^/ tpcj) + 1 


4e, 2 

((A - ^i)^p + (7 - 6)(j)) + C(a - /3) 


o{rf), 


\iva{u-£j)y/ip-(j)-+ ly/'il)+(j)++ ld{u)W ^^ = 7 EW+ 1 )-C^)+ 


{u - Vi){u + Vi) 


Vl 

2 


- {£%i^(t> +1) - ^2)^ + \{£]{.i’4> +1) - ^^)E yE" 

i=i S' s ^ S' ^k 


- y(V"/>+ 1) - ^ + ^((A + /7)S+ (7 + S</') - (a +/3)ejVEn+ 




2y/i)(i) +1 


E , fl 

4£, ^ 2 

((A - ^i)tp + (7 - S)(j)) - ^{a - (3) 


+ 0 (^ 2 ). 
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Finally, the sum of the first two terms in ( 28 ) can be expressed as 


= V 


lim {u — ej)\/^p~(j>~ + l\/'tp+(j)+ + 1 

'j 


c c. 

, .YJ {u-Vi + r]){u + Vi+T]) ■,Yl{u-Vi-rii)(u + v,-'q) 






(u - Vi){u + Vi) 


^ £2 


E 


_^ 26 '- 3 

~ E 2 '^2+7 — ej{'lp(j) +!)—(« + + 1 + 


2 2 , , ^ 

,k^j J ^ 1—1 J 


iSj 




+ 0 ( 77 ^). 


2 -v/V '0 + 1 

The third term of ( 28 ), reproduced here for convenience. 


lim {u — Sj ) 


n —:-7 / , , —- , 2 2 /^\~'' 7 ^/ 4 )(('w — £i)^ — 77^/4) 

+ 1^^+^+ + 1 C 77^ - 77 V 4 n - (2 2 2 - i - — 

M [U^ - £f)(U^ - vf) 


is computed as follows. First, expand the product in powers of 77: 

A ((u + Sj)^ - 77 ^/ 4 )((u - - 77V4) _ (u + £if(u- Sif - (772/4) ((u + £i )2 + (u - £if) + 0(772) 

11 (y2 _ ^2^^y2 _ .y2^ 11 


,,2 ^2 „2 7: 2 p 2 1 

n 77 — gj _ T]_ U — Si i 

' 7/2 / 1 . 2 ^ r ^ 


Z =1 


i,k\i^k 


2 = 1 


— vf u‘^ — vi \u — Ek U-\- Sk 


{u + ei){u + ei){u‘^ -vf) 
0 ( 772 ). 


Thus, 


lim-£,)[] 


((■u + £i)2 - r72/4)(('U - £j)2 - 772/4) 772 .^£2- £2 2£,- 


(^^2-£2)(772-u2) 


^ _:L Viz 
4 _ , 


4 2 :^ £ V vf £2 - 7;2 

i^j J * 7 J 


7 ( 77 ^). 


This term already gives the multiple of if', so we just need to consider the constant contribution from the other 
multiples. Consider 

\/+ 1VtZ+V + 1 cl = ( 7 /“(/+ + (/“ 7 /)+) + 2 - 2yVVlMVE^^^I^|r,=0 = 


Itj—O 


= 27//) + 2 - 2(7/>(/+ 1) =0 

Thus, there will be no contribution in the eigenvalues from the third term in ( 28 ). 
Finally, we obtain the eigenvalues of the conserved operators ( 24 ) as 


_(£ 2 (V,</+ 1 )-^ 2 ) 

■^3 ~ 


^ £2 
^3 


£,- 


^£3 


2 v1 V+T 


E 

((A - 77)7/+ (7 - (J)/!). 


2£2 


_ p2 _ 2 “ E -2 J 2 4 “ +!)-(« + p)£j 1 + 

k^j ^3 S' J ( 29 ) 


Remark 4.1. In view of ( 26 ), by setting j 3 = 'f! = 4 > = 5 = yL = ( = Q in ( 29 ) we deduce the eigenvalues of t* to 
be 

C 2 C r, 2 

po) 


iW, £ t S £ - « 


4.3 Bethe Ansatz Equations 

The eigenvalue expression for A(7i) given in ( 28 ) is undefined for u = Vk, for each k = 1 , 2 ,.. .£. Assuming that 
the Vk are all distinct, analyticity of A(77) requires that lim A(77) must be finite for each k = 1 , 2 , ...,£. This 

U—>-Vk 
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requirement equates to evaluating the residue of A{u) at u = and the resulting constraints on the Vk are 
referred to as the Bethe Ansatz Equations (BAE). The BAE are equivalent to 

lim {u — Vk)A{u) = 0, (31) 

U—¥Vk 

Compute (31) from (28): 

^ + 1 + + l) + 1 + ?"'’ + I+ l) X 

c ^ c 

X TT 7 -^—7^777- ^^—ThT “ Vi + T]){vk + + Tj)- 

(Vk -si+ v/^){vk + ei+ r]/2) 

- ^ (vk\/i’~4>~ +1 - - I+1) (vk\/i’~''4>~'' +1 - C'*' - +1^ X (32) 

c ^ c 

X n 7-77,77-^-77,7 + 

i=i w 2)(?7fe + £/ - 77/2) 

+ ((7/>"<)i+ + <?i"V'’^) + 2 - 2y/'ilj-(j)- + l^/i;+(l>+ + l) = 0. 

Remark 4.2. One can also compute the BAE from 

lim (n + nfc)A(M) = 0 . 

U —>- — Vk 


This gives the same expression (32). 


4.4 The quasi-classical limit of the BAE 

Let us expand the BAE (32) in the powers of p. Start with 
£ . £ 


f=i + Tll‘i){vk + £/ + 77/2) 


n 

c 

n 


1 


n 

1=1 

£ 

n 


1 


2 ^ 


i=i 

£ 


1 


1 


[vk -£i- r]/2){vk + £/ - 77/2) v% - Sj 


' + 7 E 


Vk - £i Vk+ £i 

1 1 

+ 


Vk -£i Vk+ £i 


0(77), 

■0(77), 


Thus, 


'[[{vk-Vi+p){vk + v^ + p) =Y[{vl-v^) I 1+77X1 
i^k \ i^k 

C C / C 

Yl{vk-Vi-p){vk + v,-p)=Yl{vl-v^) I-77X 


i^k 




i^k 


1 1 

-^- 

Vk -Vi Vk + Vi 

1 1 

-^- 

Vk -Vi Vk + Vi 


n 


y (». -.,+,/ 2 )(„,+.,+,/ 2 ) n(»‘ - ». + 7)(« + .. + ^) = 


£ ^ £ 


=n7+72n+-+i‘+7|E7^+-E7+ 


Vk 


0(77), 

0(77). 


/ — I ii^k 


j^k ^1 m=l ^ra 


0 ( 77 ), 


n 

c 

n 


£ £ 

117 - 7T\( -^- j^\\{vk - Vi - p){vk + V, - p) = 

i=i “ 77 / 2 ) (£fe + £/ - 77 / 2 ) 


£ ^ £ 


/ t: 


n(^fc-+") 


2 2 

'I’r — £7 

/=1 ^ ^ i^k 


Vk 




_ p 2 

1=1 ^k 


0(77). 
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One can check that the first order contribution (i.e. first order in powers of rj) from the third term in the sum on 
the left hand side of (32) is zero: 


f)+ + + 2 - 2 a/V'"</>" + + l) = 


= + ri{X-tp + 5(1)) + +-f(j))) + 2 - 2(V’<() + 1) - 2(V^^ ^ 


The first order contribution from the other two terms in (32) also gives zero: 

c 

(ffcW’<?^ + l-C) (vk\/H)Tl + n „ 2^^2 


2ri 

Vk 


2r] 

Vk 


c ^ c 

~2 2 

VT — ST 

1—1 ^ ^ i^k 

c ^ c 




Vt — ST 

1—1 ^ ^ i^k 


Thus we have to expand the BAE up to the second order. Start with 

^,p 4, +l=V^^)l + 5--8 W + 1)^ 

1 rnrrx-f) , 1 ?7(AV’+ 7<(') + »7^A7 lr]'^{Xi) + 'y(j))^\ 2^ 

+l='/^^(' + 2-vTTl-8 »4,+ lf ' 

Using this, we now calculate the contribution from the third term of (32): 


SX + — 2('0^ + 1) 

{{X - n)^) + {-f - 6)(j)) 


+ X- 


1 (liip + d(l)){X'tp + j(j)) ^ 1 /iJ 1 (fj,)/) + 5(5)'^ ^ 1 A 7 1 (A-)/) + 7 ^) 
4 (^(f) + 1)2 

-{'^-5){X- ^i). 


21 


2 i/jc)) + 1 8 (■;/'(/)+1)2 2 ■0^ + 1 8 (00+1)2 


4(00 + 1) 

Using similar techniques, the first term of (32) gives the rf contribution 


/=1 ^ ^ i^k 


(vk v/00 + 1 - C) 


Ufc A0 + 70 000 + 1 

2 000 + 1 2 


+ (ufc \/00 + 1 + ^ 


Ufc ^0 + 00 
2 000+1 


000+1 

2 




K( 00 +i)-^ 2 ) I E 






and the second term gives the rj^ contribution 

(vk 000+1 + C 

+ (vk\/tp(j)+ 1 -■?) 


n 0 - 0 

1—1 ^ ^ i^k 


2 000 + 1 

Ufc ^0 + 00 _ o _ 0 V’ 0 + 1 
2 000 + 1 ^ 2 


Vk A0 + 70 000+1 

— a--- 


/ £ 


-(7w0+i)-niE003-E0^ 






Summing up all terms we obtain 
c , c 


2n;0T5n+-+ 2 (« + fflvWTi + 2W + i)-««1^0|±ip+> + 


<11-^ _ c*- 




IE ^ 

E Ml 

\h 

( ^ 'iP' _ p 2 1 

m -1 / _ 


000 + 1 

((A-/x) 0+(7-j)0)' 
4(00 + 1) 


- (7 - 5){X - ^) = 0. 
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Thus, we obtain the following BAE in the quasi-classical limit (keeping in mind the parameters (j), ijj, a, /3, 7 , S, 
A and fj, are defined in (20)): 


+ (vli-tpcj) + 1) - ;; 2^=-2 

\i:^k ^ 




= j i (7-<5)(A-m) - 

4 


((A-m)V' + (7-^)<(') ) ULii^k-sf) 


4 ( 7 / 1 ^ +1) 




(33) 


By setting P = ip 
satisfy the BAE 


= (5 = /r = ^ = 0in (33) we deduce the Bethe roots {vk : k = 1, ■■■,£} appearing in (30) 


G; A 1 



E 



7A UtM-ef) 

uU(^i-vfy 


(34) 


Remark 4.3. Previous studies have used a correspondence between BAE and differential equations, through a 
generalised Heine-Stieltjes problem, as a route to numerically solve BAE for a wide range of models [17,26-31]. It 
is interesting to observe that the BAE (3f) have an interpretation as an inhomogeneous, generalised Heine-Stieltjes 
problem. Define the polynomials 


c 

Q{x) = Ylix-Vj), 

i=i 

Pix) = Ylix-ej). 

i=i 


The BAE (34) are equivalent to the condition 

vlP{vl)Q”{vl) + ((a + l)P(vl) - vlP'{vl)) Q'{vl) = ^ [PpvDf . 

It follows that Q{x) satisfies an inhomogeneous, linear, second-order differential equation 

xP{x)Q''(x) + ((a + 1)-P(a;) — xP'{x)) Q'{x) + V{x)Q{x) = ^ [R(a;))]^ 

where V(x) is a Van Vleck polynomial of order C. A .similar correspondence also applies at the level of the BAE 
(33). 


4.5 Eigenvalues of the Hamiltonian 

Recall A*, given in (30), is the eigenvalue of the conserved operator t* given in (26). To compute the eigenvalue 
of the Hamiltonian (27) consider 


c c 


Es'E*= E a^-EE 


j=i 


^— 2 - q: g,- ^ = 

_ J J < J J 

j,k:k^j i k j i 


= -2 —2 ~ q^e^ 

-E -e] ^ 

i,]=i J j=i 


From the BAE (34) we find 


E 1 ^ G + 1 E 2 _ 7 A - e|) 

2-^ ^2 — ^2 2-^ 7,2 — 


VT ^.vt-vl 4u2nE(^2_^2)' 


k^i 
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Thus, 


2 


4 


This leads to 


£ £ 


2 (a + l)^z;j + ^ ^2 _ 2 2 rr-^ r 2 2', 

i=l J = 1 i=l i,k-.i^k^^ 1=1 ^fe) 


i=l 

£ 


= 2(a + l) Vu-2- 

'^--2\*_r,/ I i'i'^,.-2 7^1 nj=iy“£i) y 

\ 2 (a+i)yv^ y„,2T-r£ ,9 „,2i ^y* 


.-2 


i=i 


2 


i=i 


Implementing the change of variables Zj = £ ■ ^, yi = and setting A = —7 we obtain the expression 


z. 11, 1^2, r^y n,=i(i-2/i^, ) 

£; = (i+G)yy.--yzy-y-^-—^ 

1=1 ^ k=i ^ i=i UkM^-y^yk ) 


(35) 


for the eigenvalues of the Hamiltonian (3) subject to the BAE obtained from (34): 


1+G-i+y ^^+y = -^2 

^ y^-yk i^yk- zf G^yk Y\f,^{i - 


nya-y^y) 

„-11 


i^k 


Z =1 


Vky, 


(36) 


5 Conclusion 

We have shown that the Hamiltonian (1), describing a p +ip pairing model interacting with its environment, is an 
integrable model. By mapping the Hamiltonian to the spin operator formalism (3), we found through a derivation 
of the model by the Boundary Quantum Inverse Scattering Method that the energies of (3) are given by (35) 
subject to the Bethe Ansatz Equations (36), and the operators (26) are conserved with eigenvalues given by (30). 
It is anticipated that this exact result will allow for a detailed analysis of the model in future studies. 

There are two examples in the literature [32, 33] of related pairing models interacting with a single bosonic 
degree of freedom, where the boson-fermion interaction has the form 

^ X] ~ iky)c-]^ci^b^^ . 

k 

This is in analogy with the system-environment interaction incorporated into (1), however as mentioned in the 
Introduction for these models density matrices of the pairing model generically exhibit entanglement with the 
bosonic degree of freedom. These two models have vastly different ground-state behaviour, with [33] exhibiting 
features which are qualitatively similar to those of the p + ip pairing Hamiltonian while those of [32] are not. It 
will be interesting to see what (3) has to offer in this regard, and what the consequences are of not accommodating 
entanglement with the environment. Another important line of future research is to extend the body of results on 
exact form factors and correlation functions for the p + ip pairing model [12,16], and the same model interacting 
with a bosonic mode [32,34], to the integrable generalisation (1). 
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